Influence of van der Waals Forces on the Adsorption Structure of Benzene on Silicon 
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Two different adsorption configurations of benzene on the Si(001)-(2x 1) surface, the tight-bridge 
and butterfly structures, were studied using density functional theory. Several exchange and correla- 
tion functionals were used, including the recently developed vdW-DF functional, which accounts for 
the effect of van der Waals forces. In contrast to the PBE, revPBE and other GGA functionals, the 
vdW-DF functional finds that, for most coverages, the adsorption energy of the butterfly structure 
is greater than that of the tight-bridge structure. 
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In the quest for reduced-sized transistor chips the com- 
bination of organic molecules with silicon-based technol- 
ogy is of increasing importance. The ability to manipu- 
late organic molecules on surfaces is developing rapidly 
and an understanding of the structural and transport 
properties of adsorbed molecules is essential [H, 0, 0, 0, 0] . 

Silicon is not only important for technology but also 
demonstrates the versatility of the covalent bond. Cova- 
lency makes bulk silicon, diamond, and graphene layers 
very strong. At the same time, it can produce a multi- 
tude of competing atomic structures when spatial restric- 
tions are released. For instance, the Si(lll) and Si(OOl) 
surfaces show quite different properties, including radi- 
cally different types of reconstructions. While the cova- 
lent bond is typically very strong, the energy differences 
between such reconstructions can be small [6| . Their rel- 
ative stabilities and the influence from adsorbates are in- 
teresting issues to understand. In particular, there may 
be situations where the weak van der Waals (vdW) force 
can significantly influence the structures. The benzene 
molecule demonstrates several types of internal bonds, 
typically interacts with other molecules via vdW forces, 
and is an important model unit for several classes of large 
molecules like DNA Q ■ The adsorption of benzene on Si 
is obviously very interesting, both because of the versa- 
tility in bond types and the wide ramifications, including 
technological ones, such as molecular electronics. 

Adsorption of benzene on the Si(001)-(2x 1) surface 
has been studied extensively, however, without unani- 
mous results. There is agreement on the two most stable 
structures, shown in Fig. [U but to date there is no agree- 
ment which is the stable one. The butterfly (BF) struc- 
ture is adsorbed on top of a single dimer and has two 
symmetry planes: along and perpendicular to the dimer 
row. The tight-bridge (TB) structure is adsorbed across 
two dimers and has one symmetry plane along the dimer 
row. 

To differentiate between the adsorption structures and 
their bonding and symmetry properties, a variety of ex- 
perimental tools have been employed. Thermal desorp- 
tion and angle-resolved photoelectron spectroscopy were 
used to investigate the electronic structure and symme- 
try of benzene on Si [|[ and a single dimer structure was 







FIG. 1: (Color online) Butterfly (BF) (top), and tight-bridge 
(TB) (bottom) structures of C 6 H 6 on Si(001)-(2x 1). LHS: 
viewed along [110]. RHS: viewed along [110] dimer rows. 



observed, supporting the BF configuration. Near-edge 
X-ray-absorption fine-structure found the benzene to be 
symmetric with respect to the dimer axis, ruling out the 
TB structure 0] . These findings were supported by opti- 
cal spectroscopy data [l(J , which found that benzene ad- 
sorbs on top of a single dimer rather than on the bridge 
site between two dimers. However, according to scanning 
tunnelling microscopy (STM) studies [HEl, thc bcn " 
zene molecule adsorbs initially in the BF structure but 
this is observed to be metastable with respect to a bridg- 
ing configuration. With the STM tip one benzene struc- 
ture can be converted to another, with the conversion- 
energy barrier estimated to be 0.95 eV. A high-resolution 
study [IH suggested that the adsorption geometry de- 
pends on coverage [37| and showed a bridging structure 
is favored at low coverages while at high coverages a sin- 
gle dimer structure (BF) is more stable. This would ex- 
plain the discrepancy between the experimental results 
of Refs. 0, 0, [TO that supported the BF structure, as 
they were performed at saturation coverage, while the 
STM experiments were carried out at a low coverage of 
0.044 ML. 

Adsorption energies, i? a ds, have been calculated using 
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standard densit y fu nctional theory (DFT) and all of these 
studies [U, [llllil support the TB structure. An MP 2 
(second order perturbation/quantum mechanics molecu- 
lar dynamics) cluster calculation [ItJ is the only calcu- 
lation which supports the BF structure (.E a ds=L04 eV). 
Unlike standard DFT calculations, MP2 methods include 
vdW forces but because they are computationally heavy 
only small systems can be treated with this method. Si 
(001) is an extended material and it is not clear whether 
small clusters can give accurate adsorption energies. Fur- 
thermore, the cluster geometry corresponds to a low cov- 
erage situation so cannot be compared with saturation 
coverage results. 

The vdW-DF functional [l8| was developed to account 
for the effect of vdW forces in DFT. It has been shown to 
give accurate results for molecular systems, such as ben- 
zene dimers [lE and recently it has been success- 
fully applied to systems with covalent bonding present 
pH |22| . [23} . It is currently believed that vdW forces are 
only important in physisorbcd systems but in this pa- 
per we demonstrate that this is not the case. By using 
vdW-DF in DFT calculations we show that in most of 
the studied configurations vdW forces stabilise the BF 
structure, which is surprising since the TB structure has 
more covalent bonds and would be expected to be more 
strongly bound to the surface. 

Standard DFT calculations were performed with Da- 
capo using the PBE form [2E HE H3 of the 
generalised-gradient approximation (GGA) and ultrasoft 
(US) pseudopotcntials. The plane- wave energy cutoff was 
400 cV and the Brillouin zone mesh used was equivalent 
to 4x4x1 Monkhorst-Pack fc-point sampling for 0.5 ML 
coverages. We used 9-atomic layer Si slabs, «15 A of 
vacuum and a lattice constant of 5.47 A. To save compu- 
tational time adsorption on only one side of the slab was 
considered. The bottom layer of the slab was fixed in the 
bulk Si positions and passivated with two H atoms per Si 
atom. The slabs were electrostatically decoupled along 
the c-direction. All relaxations were considered complete 
when the forces were less than 10 meVA -1 . 

Chemisorption energies are known to be sensitive to 
the choice of the exchange and correlation (xc) func- 
tional. It has been shown that for chemisorbed atoms and 
molecules on transition metal surfaces the revPBE func- 
tional [28j] generally gives better results than the PBE 
functional, which tends to overbind the molecules [2^|. A 
comparison between the present E^ds and those of pre- 
vious studies is shown in Table [T] The PBE and PW91 
[30II3H xc functionals give similar E a ds, as expected, but 
are around two times higher than the revPBE E^ds- 

The additional binding energy due to vdW forces is 
calculated using the post-GGA total en ergy vdW density 
functional (vdW-DF) described in Ref. [l8|. To minimize 
any artificial exchange binding which can be mistaken 
from vdW-binding, the vdW-DF employs the revPBE 
form for the exchange description [33|, [34|] . The vdW-DF 
replaces the GGA correlation and divides the correlation 
into a shorter ranged and a longer ranged part. The 
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Pseudopotential GGA 


Reference 


0.66 0.47 


US 


revPBE 


Present 


1.16 0.89 


US 


PBE 


Present 


1.24 0.99 


us 


PW91 


Present 


1.25 1.00 


PAW 


PW91 


[32] 


1.05 0.82 


US/NC 


PBE 


U5] 



TABLE I: _E a d s of benzene in the TB and BF geometries for a 
coverage of 0.5 ML. The pseudopotentials (psp) and xc func- 
tionals used are also shown. 



first part is approximated by the LDA, while the lat- 
ter (E^) includes the important dispersive interactions. 
E^ is nonlocal by construction and, consistent with the 
approximation for the shorter ranged correlation, it is 
constructed to vanish for a homogeneous system. In the 
vdW-DF the total energy reads: 

^nvdW-DF _ ^.revPBE _ ^irevPBE _|_ ^LDA _|_ j^c 

with all terms obtained on the basis of self-consistent 
semi-local PBE DFT calculations. The nonlocal correla- 
tion can be written in the simple form [l8j : 

Ki = lJ dr J dr'n(r)0(r,r>(r'), (2) 

with a density-density interaction kernel, (f>(r, r') derived 
from the many-body response of the weakly inhomogc- 
neous electron gas. The kernel </>(r, r') can be tabulated 
in terms of two dimensionless parameters, a scaled sep- 
aration D = |r — r'| (go + <7o)/2 and an asymmetry pa- 
rameter S = (go — Qo)/{<lo + <?o)j w here qo(r) is a local 
parameter that depends on the electron density and its 
gradient at r. More information on the explicit form and 
deriviation of the kernel can be found in Ref. 

E^i is to be evaluated for a periodic system and we use 
the scheme developed in [IE [35|] to evaluate the integral. 
In short, we let r in Eq. ([2]) run through all electron den- 
sity points within the DFT unit cell. The primed spatial 
coordinate on the other hand, includes the electron den- 
sity from the surrounding lateral repeated images until 
convergence is obtained. We have found the integral to 
be converged to include in total 5 (3) unit cells in each 
direction for the largest (lowest) coverages. 

The post-GGA version of the vdW-DF functional im- 
plemented in this article does not allow for any additional 
geometric and electronic relaxation beyond the PBE de- 
termined adsorption structure. Thonhauser et al [2l| im- 
plemented the vdW-DF functional self-consistcntly and 
for noble gas and cytosine dimers the difference between 
the self-consistent and non-selfconsistent energies was 
negligible. This study also looked at the effect of the 
vdW-DF functional on bulk Si. In this expected, 
the contribution of the E^ energy was negligible and the 
lattice constant and total energy were similar to those ob- 
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taincd using PBE. This indicates that our results, which 
use the non-selfconsistent implementation, are reliable. 

Due to numerical convergence issues of the nonlocal 
integral 0, HU, special care is taken when correcting 
for the nonlocal energy in the adsorption system, lead- 
ing to additional steps in the evaluation procedure of 
vdW-DF. First, £? a d s between the benzene layer and the 
surface (BLS) is evaluated, and second, the intralayer 
molecular-molecular (IMM) energy associated by isolat- 
ing the benzene from the benzene layer is calculated. To 
ensure maximal error cancelation in the evaluation of the 
above intermediates, the isolated layer, the silicon surface 
and the benzene molecule are fixed in the PBE adsor- 
bant atomic configurations and the spatial separation of 
the density fast-fourier transform grid is kept constant 
in all calculations. Finally, the contribution of unfold- 
ing (UF) the isolated molecule and the surface to their 
PBE relaxed structures is calculated. We choose in ac- 
cord with Ref. [HI to describe this final energetic con- 
tribution within the PBE approximation. This is jus- 
tified, as all structural information has been obtained 
self-consistently within the PBE functional, and further- 
more, PBE is, in contrary to vdW-DF, designed with the 
energetics of the internal atomic binding in mind. 

The detailed data for the 0.5 ML case are given in 
Table [III with the vdW-DF total energy separated into 
the nonlocal correlation including dispersion forces 
and the remaining ^ vdwo part. We will briefly discuss 



case is about 0.24 eV more attractive. Thus, compared to 
revPBE correlation energy the BLS part of the vdW-DF 
correlation increases the BF E a d s by 0.26 eV more than 
in the TB case. This is in accord with the expectation 
that the BF structure has a larger vdW-like interaction 
than the more covalently bound TB structure. 

IMM interactions: We find E^ to be small and attrac- 
tive (« 0.05 eV) while E vdwo is slightly repulsive in both 
the TB and the BF configurations. The IMM E vdwo is 
found to be slightly more repulsive for the TB configura- 
tion than for the BF, such that the net vdW-DF energy 
is slightly repulsive for the TB case and weakly attrac- 
tive for the BF case, closely resembling the interactions 
obtained in the semi-local revPBE and PBE approxima- 
tions. As the coverage is reduced, the vdW-DF IMM in- 
teractions are even less pronounced and can be described 
with the semi-local DFT functionals to within 0.02 eV. 

Overall, the change in i? a( j s calculated with vdW-DF is 
substantial, in particular for the BF structure. The main 
contribution to £ a d s comes from between the surface 
and the benzene layer. The attraction between the ben- 
zene molecules in the layer is minimal. This behaviour is 
similar to that of a cytosine dimer [2lj | where the repul- 
sive E vdwo term is compensated by the attractive Ef^ to 
give an overall binding energy of around 0.3 eV. 

We will now discuss the effect of coverage on the ad- 
sorption energy. The adsorption energies for various cov- 
erages are shown in Table [TTT1 Both PBE and revPBE fa- 



TABLE II: vdW-DF -E a d s and its contributions, are shown for 
a coverage of 0.5 ML. The standard DFT PW91, PBE and 
revPBE results are shown for comparison. All energies are in 
eV. 





PW91 PBE revPBE vdW-DF 


TTlVdWU 7T1C 
& All 


BLS (BF) 
BLS (TB) 


3.51 3.45 3.03 3.32 
8.51 8.45 7.98 8.02 


1.96 1.36 
6.56 1.46 


IMM (BF) 
IMM (TB) 


0.00 -0.03 -0.05 0.02 
-0.03 -0.06 -0.07 -0.02 


-0.04 0.06 
-0.07 0.05 


UF (BF) 
UF (TB) 


-2.52 -2.52 -2.51 (-2.52) 
-7.23 -7.23 -7.24 (-7.23) 




£ads (BF) 

-E a ds (TB) 


0.99 0.89 0.47 0.82 
1.24 1.16 0.66 0.77 





both of the vdW-DF contributions that contain nonlocal 
components (BLS and IMM) for the 0.5 ML case, and 
compare these to the corresponding PBE and revPBE 
results. 

BLS interactions: E^ for the TB case is found to be 
0.1 eV more attractive than in the BF case. This stems 
from the fact that the TB structure is more closely bound 
to the Si surface, and accordingly we find «0.1 eV dif- 
ference for all coverages. However, E^ must be added to 
^vdwo tQ find the v dW-DF energy. The E vdwo of the 
TB configuration is at first glance considerably larger in 
the BF case. However, comparing the combined vdW-DF 
energy with the corresponding revPBE value, we see that 
the TB case has almost the same energy, while the BF 



TABLE III: Variation of adsorption energy with coverage for 
the BF and TB structures. The PW91, PBE and revPBE 
results are also shown for comparison. All energies are in eV. 





PW91 


PBE revPBE 


vdW-DF 


BF-0.5 


0.99 


0.89 


0.47 


0.82 


BF-0.25a 


1.02 


0.93 


0.51 


0.82 


BF-0.25b 


1.02 


0.93 


0.61 


0.84 


BF-0.125 


1.04 


0.96 


0.54 


0.84 


TB-0.5 


1.24 


1.16 


0.66 


0.77 


TB-0.25a 


1.23 


1.16 


0.67 


0.74 


TB-0.25b 


1.33 


1.25 


0.76 


0.86 


TB-0.125 


1.31 


1.24 


0.75 


0.82 



vor the TB configuration and, as expected, PBE predicts 
the larger adsorption energy. On the contrary, vdW-DF, 
with its account of the dispersion interactions, predicts 
the BF configuration to be slightly favored. For a cov- 
erage of 0.25 ML two supcrccll orientations arc possible 
and are denoted (a) and (b) with primitive lattice vec- 
tors (220)(110)(006) and (110)(220)(006), respectively. 
In general only a minimal coverage dependence is found 
as changes in E^ are almost cancelled by the correspond- 
ing changes in E vdwo . Correspondingly, the vdW-DF 
adsorption energies in the BF case only show a small in- 
crease with increasing coverage. For the TB structure, 
the coverage dependence is almost constant except for 
0.25b coverage, which has a pronounced preference. This 
effect is also observed for the semi-local functionals. For 
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all coverages the semi-local functional favor the TB con- 
figuration. In contrast, vdW-DF predicts that the BF 
structure is stable for all cases, except 0.25b where the 
TB is lower in energy by only 0.02 eV. 

In the 0.125 ML case the nonlocal IMM interactions are 
less than 0.01 eV so the monomers can be regarded as iso- 
lated. The difference between the BF and TB configura- 
tions is thus 0.02 eV in the isolated case, which is slightly 
less than the 0.08 eV difference found in a correspond- 
ing MP2 calculation (l7| . To compare with experimental 
data we used the Redhead equation (3(| to estimate the 
adsorption energies based on thermal desorption spectra 
in Ref. Q. For a peak temperature of 432 K, a heating 
rate of 5 Ks _1 and a pre-exponential frequency factor 
between 10 12 -10 16 s _1 , we obtain an adsorption energy 
in the range 1.06-1.40 eV. The difference between exper- 
iment and vdW-DF energies is due to the uncertainty 
in the GGA exchange energies [p| , which can be seen in 
Table HH The use of PBE with the vdW-DF increases the 
adsorption energies to within the experimental range but 
the uncertainty in the exchange energy is large enough to 



mask the small energy differences between the two struc- 
tures, particularly for the lower coverages. 

In summary, we have demonstrated that standard 
DFT adsorption energies are significantly dependent on 
the xc functional and, furthermore, the inclusion of vdW 
forces makes a qualitative difference to the results. Stan- 
dard DFT with PBE and revPBE functionals finds that 
the TB structure is always stable, whereas vdW-DF DFT 
calculations find that, for some coverages, vdW forces 
stabilise the BF structure. These results have significant 
implications for many DFT studies as vdW forces are 
generally considered to have a negligible effect on cova- 
lently bonded systems and are usually ignored. 
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[37] In this paper one monolayer (ML) corresponds to one 
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